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'We investigate the formation of the double pulsar PSR J0737-3039 and examine its most likely 
progenitors, taking into account the most recent and all currently available observational constraints. 
We show that the most likely kick velocity and progenitor parameters depend strongly on the 
consideration of the full five-dimensional probability distribution function for the magnitude and 
direction of the kick velocity imparted to pulsar B at birth, the mass of pulsar B's pre-supernova 
helium star progenitor, and the pre-supernova orbital separation rather than marginalized one- or 
two-dimensional distributions for the kick velocity and progenitor mass. The priors that enter the 
analysis are the age of the system, the minimum helium star mass required to form a neutron star, 
the transverse systemic velocity, and the treatment of the unknown radial velocity. Since the latter 
cannot be measured observationally, we adopt a statistical approach and use theoretical radial- 
velocity distributions obtained from population synthesis calculations for coalescing double neutron 
stars. 'We find that when the minimum pre-supernova helium star mass required for neutron star 
formation is assumed to be 2.1 Mq, the most likely kick velocity ranges from 70 kms~^ to 180 kms~^. 
When, on the other hand, masses lower than 2.1 Mq are allowed as neutron star progenitors, the 
most likely kick velocity can reach very low values (as low as a few kms~^), although the majority 
of the models still yield most likely kick velocities of 50kms~^ to 170kms~^. Hence, we agree 
with Piran & Shaviv [19] that the observed system properties, including the low transverse systemic 
velocity, can indeed be compatible with low progenitor masses and low kick velocities. Equally 
important though, we show that this is not the only likely formation path of pulsar B, due to 
the role of different prior assumptions that are necessary in the analysis. Moreover, in contrast to 
earlier claims in the literature, we show that the proximity of the double pulsar to the Galactic plane 
and the small proper motion do not pose stringent constraints on the kick velocity and progenitor 
mass of pulsar B at all. Instead, the constraints imposed by the current orbital semi-major axis 
and eccentricity and the orbital dynamics of asymmetric supernova explosions turn out to be much 
more restrictive. 'We conclude that without further knowledge of the priors, the currently available 
observational constraints cannot be used to unambiguously favor a specific core-collapse and neutron 
star formation mechanism. Both electron capture and neutrino-driven supernovae therefore remain 
viable formation mechanisms for pulsar B. 

PACS numbers: 97.60.Bw,97.60.Gb,97.80.-d 



I. INTRODUCTION 

Since its discovery in 2003 [1, 2], the double pulsar 
PSR J0737-3039 has proven to be an extraordinary sys- 
tem for the study of both fundamental physics and gen- 
eral relativity [3, 4]. The system is highly relativistic 
with an orbital period of 2.4 hrs and an eccentricity of 
0.0878, making it the tightest known double neutron 
star (DNS) binary known to date. The system also har- 
bors the fastest known DNS pulsar (22.7ms), has the 
largest apsidal-motion rate of all presently known DNSs 



'Electronic address: b-wlllemsanorthwestern.edu 
t Electronic address: j-kaplan- lOnorthwestern . edu 
tElectronic address: tassosfragosanorthwestern.edu 
^Electronic address: vickyanorthwestern.edu 
^Electronic address: kbelczynanjnsu.edu 



(16.9° yr ^), and shows eclipses of the fast pulsar by 
the magnetosphere of the slow one (orbital inclination 
I = 87° ± 3°) [1, 2]. Its short merger time of ~ 85Myr 
has furthermore led to a significant upward revision for 
the merger rate of coalescing DNSs [5, 6]. 

A central topic in understanding the physics of the 
double pulsar is the study of its origin and evolution. 
According to our current knowledge (for more details see 
[7-10]), DNSs form from primordial binaries in which, 
possibly after some initial mass-transfer phase, both com- 
ponent stars have masses in excess oi ^ % — \2 Mq. Af- 
ter the primary explodes in a supernova (SN) explosion 
to form the first neutron star (NS), the binary enters a 
high-mass X-ray binary phase in which the NS accretes 
matter from the wind of its companion. The phase ends 
when the companion star evolves of the main sequence 
and engulfs the NS in its expanding envelope. The NS 
then spirals in towards the core of the companion until 
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enough orbital energy is transferred to the envelope to ex- 
pel it from the system. When the envelope is ejected, the 
binary consists of the NS and the stripped-down helium 
core of its former giant companion, orbiting each other 
in a tight orbit. If the NS is able to accrete during its 
rapid inspiral, a first "recycling" may take place during 
which it is spun up to millisecond periods. The helium 
star then evolves further until it, in turn, explodes and 
forms a NS. Depending on the mass of the helium star 
and the size of the orbit at the time of the explosion, the 
SN may be preceded by a second recycling phase when 
the helium star fills its Roche lobe and transfers mass to 
the NS [9, 11, 12]. 

Shortly after the discovery of PSR J0737-3039, Dewi & 
van den Heuvcl [13] and Willems & Kalogera [14] (here- 
after Paper I) independently derived that, right before 
the second SN explosion, the binary was so tight that 
the helium star must have been overflowing its critical 
Roche lobe. This conclusion, for the first time, strongly 
confirmed the above formation channel for DNS binaries. 
The observed 22.7 ms pulsar (hereafter PSR J0737-3039A 
or pulsar A) then corresponds to the first-born NS, and 
its 2.8 s pulsar companion (hereafter PSR J0737-3039B or 
pulsar B) to the second-born NS. 

In a follow-up study to Paper I, Willems et al. [15] 
(hereafter Paper II) extended the stellar and binary evo- 
lution constraints on the formation of pulsar B to in- 
clude constraints based on the system's position and mo- 
tion in the Galaxy. While the addition of proper motion 
measurements to the list of observational constraints, at 
first sight, added a significant piece of information on 
the formation and evolution of the double pulsar system, 
the unknown direction of the proper motion and the un- 
known radial velocity significantly hindered the deriva- 
tion of more stricter progenitor constraints than obtained 
in Paper I. 

Since Paper II, the available observational constraints 
have undergone significant revisions. Most notably, 
proper motion estimates, which form a vital part in con- 
straining the kick imparted to pulsar B at birth, have 
decreased from ~ 141kms~^ down to ~ 30kms~^ or 
less [4, 16, 17]. In addition, rough limits on the spin- 
orbit misalignment of pulsar A have recently been de- 
rived by Manchester et al. [18], which again add to the 
constrains on pulsar B's natal kick velocity. Prom a theo- 
retical point of view, it has furthermore been argued that 
pulsar B may have had an unusually low- mass progenitor 
(--v 1.45 Mq) which formed a NS through a new type of 
gravitational collapse accompanied by an unusually low 
kick velocity (< SOkms^i) [19, 20]. 

Our aim in this paper is to investigate the system's 
recent evolutionary history and examine the necessity of 
this new type of gravitational collapse to explain the for- 
mation of pulsar B. We incorporate the most up-to-date 
observational information on PSR J0737-3039, including 
the 30kms~^ upper limit on the proper motion derived 
by Kramer et al. [4]. A significant new clement com- 
pared to Papers I and II, which focused on the extent of 



the available parameter space, is the examination of the 
full multi-dimensional probability distribution functions 
(PDF) for the magnitude and direction of pulsar B's na- 
tal kick velocity, the mass of pulsar B's pre-SN helium 
star progenitor, and the orbital separation right before 
the formation of pulsar B. Other significant improve- 
ments are the incorporation of theoretical distribution 
functions for the unknown radial velocity by means of 
binary population synthesis calculations for coalescing 
DNSs, the investigation of projection and marginaliza- 
tion effects of the multi-dimensional PDF, and the ex- 
amination of the role of the prior assumptions. 

The plan of the paper is as follows. In § II, we briefiy 
summarize the observationally inferred properties of the 
double pulsar relevant to this investigation. In § III. we 
present the basic assumptions and outline the method 
used to derive the formation history of pulsar B. Details 
of the various steps adopted in the investigation and their 
overall role in the analysis are presented in § IV and § V. 
In § VI, we present the most likely progenitor and forma- 
tion parameters for PSR J0737-3039B and systematically 
explore the sensitivity of the parameters to the adopted 
assumptions. The final two sections, § VII and § VIII, 
are devoted to a comparison with previous investigations 
and summarizing remarks. 



II. OBSERVATIONAL CONSTRAINTS 

The main properties of the double pulsar system and 
its component stars relevant to this investigation are 
summarized in Table I. The principal differences with the 
observational constraints used in Paper II are the consid- 
erably lower estimate for the systemic velocity perpendic- 
ular to the line-of-sight, and the availability of limits on 
the spin-orbit misalignment angle of pulsar A. 

The kinematic properties of PSR J0737-3039 have un- 
dergone significant revision since Ransom et al. [16] used 
interstellar scintillation measurements and inferred a ve- 
locity component perpendicular to the line-of-sight of 
~ 141 kms~^. Coles et al. [17] derived a reduced velocity 
of ~ 66kms~^ by incorporating anisotropics of the in- 
terstellar medium in the interstellar scintillation model. 
This reduced velocity, however, strongly depends on the 
adopted anisotropy model. Kramer et al. [4], on the other 
hand, used pulsar timing measurements to derive a firm 
model-independent upper limit of 30 km s^^ on the trans- 
verse velocity. In view of these revisions, we update the 
kinematic properties of PSR J0737-3039 derived in Pa- 
per II, adopting the most recent constraint of 30kms~^ 
as an upper limit to the velocity component transverse 
to the line-of-sight. 

The kick imparted to pulsar B at birth is expected to 
tilt the orbital plane and misalign pulsar A's spin axis 
with respect to the post-SN orbital angular momentum 
axis [e.g. 22]. The degree of misalignment depends on 
both the magnitude and the direction of the kick, and 
therefore yields a valuable piece of information on pul- 
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TABLE I: Physical properties of PSR J0737-3039. 



Parameter 


Notation 


Value 


Reference 


Distance (pc) 


d 


600 


[1,2] 


Galactic longitude (J2000) (deg) 


I 


245.2 


[1,2] 


Galactic latitude (J2000) (deg) 
Proper motion (kms~^) 


b 


-4.5 


[1,2] 


Vt 


< 30 


[4] 


Spin period of pulsar A (ms) 


Pa 


22.7 


[1,2] 


Spin period of pulsar B (s) 


Pb 


2.8 


[1,2] 


Mass of pulsar A (M©) 


Ma 


1.34 


[1,2] 


Mass of pulsar B (Mq) 


Mb 


1.25 


[1,2] 


f:nFtT'FtP^"AT*1Qi"ip P (TP ("IT T^lll'JFtT' A 1 A/r'\7'"ri 

Xw'iicii cK^uci loijiij *-''o jJU-iodi -rt- I ivi_yi ( 


TA 


210 


fl 21 


Characteristic ago of pulsar B (Myr) 


TB 


50 


[1,2] 


Spin-down age of pulsar A (Myr) 


Tad 


100 


[1,2] 


Spin-orbit misalignment of pulsar A (deg) 


A 


< 60 or > 120 


[18] 


Orbital period (hours) 




2.4 


[1,2] 


Orbital semi- major axis (-R©) 


Ac 


1.26 


[1,2] 


Orbital eccentricity 


ec 


0.0878 


[1,2] 



sar B's natal kick velocity (sec §VB for more details). 
In Paper II, we showed that for the pre- and post-SN 
orbital parameters compatible with all available observa- 
tional constraints, misalignment angles between approx- 
imately 70° and 110° are highly unlikely. Manchester et 
al. [18] derived observational constraints on the spin-orbit 
misalignment based on the stability of pulsar A's mean 
pulse profile over a time span of 3 years. The authors con- 
cluded the angle to be smaller than 60° or larger than 
~ 120°, in agreement with the theoretical predictions of 
Paper II. In the present paper, we therefore include the 
limits on pulsar A's spin-orbit misalignment inferred by 
Manchester et al. [18] in the list of available observational 
constraints. 

Among the constraints inferred from observations, the 

most uncertain parameter affecting the reconstruction of 
the system's formation and evolutionary history is its age. 
Characteristic ages, defined as half the ratio of the spin 
period to the spin-down rate, are the most commonly 
adopted age estimators, but are known to be quite unre- 
hable [e.g. 23]. In the case of PSR J0737-3039, the spin 
evolution of pulsar B is furthermore very likely affected 
by torques exerted by pulsar A on pulsar B [24], adding 
to the uncertainties of its age. Lorimcr et al. [25] there- 
fore derived an alternative age estimate by noting that, 
according to the standard DNS formation channel, the 
time expired since the end of pulsar A's spin-up phase 
should be equal to the time expired since the birth of 
pulsar B. By combining this property with a selection of 
different spin-down models, the authors derived a most 
likely age of 30-70 Myr, but were unable to firmly exclude 
younger and older ages. A third age estimate can be ob- 
tained by assuming that pulsar A was recycled to the 
maximum spin rate and calculating the time required for 
it to spin down to the currently observed value. This 
so-called spin-down age provides an upper limit to the 
age of the system of 100 Myr [1]. In view of these un- 
certainties, wc derive constraints on the formation and 
evolution of PSR J0737-3039 for three different sets of 



age assumptions: (i) < r < 100 Myr (the age range 
adopted in Papers I and II), (ii) 30 < r < 70 Myr (the 
most likely age range derived by Lorimer et al. [25]), and 
(iii) T ~ 50 Myr (the characteristic age of pulsar B). 



III. BASIC ASSUMPTIONS AND OUTLINE OF 
THE CALCULATION 

Our goal in this investigation is to reconstruct the evo- 
lutionary history of PSR J0737-3039 and constrain its 
properties at the formation time of pulsar B. More specif- 
ically, our aim is to derive a PDF for the magnitude and 
direction of the kick velocity imparted to pulsar B at 
birth, the mass of pulsar B's progenitor immediately be- 
fore it explodes in a SN explosion, and the orbital separa- 
tion of the component stars right before the SN explosion. 
Adopting the standard DNS formation channel, the pro- 
genitor of the double pulsar right before the formation of 
pulsar B, consists of the first-born NS, pulsar A, and the 
helimn star progenitor of the second-born NS. Since the 
formation of the second NS is preceded by one or more 
mass-transfer phases, tidal forces can safely be assumed 
to circularize the orbit prior to the formation of pulsar B. 
For the remainder of the paper, we refer to the times right 
before and right after the formation of pulsar B by pre- 
SN and post-SN, respectively. The subsequent analysis 
consists of four major parts. 

Firstly, the motion of the system in the Galaxy is 
traced back in time up to a maximum age of 100 Myr. 
The goal of this calculation is to derive the position and 
ccnter-of-rnass velocity of the binary right after the for- 
mation of pulsar B, and use this information to constrain 
the kick imparted to it at birth. In order to determine 
possible birth sites, our analysis is supplemented with 
the assumptions that the primordial DNS progenitor was 
born close to the Galactic plane, and that the first SN ex- 
plosion did not kick the binary too far out of it. The first 
assumption is conform with our current knowledge that 
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massive stars form and live close to the Galactic plane 
(their typical scale height is ~ 50 70 pc). The second 
assumption on the other hand neglects a small fraction 
of systems formed with high space velocities [9, 26] [52]. 
Under these assumptions, the binary is still close to the 
Galactic plane at the formation time of pulsar B. Pos- 
sible birth sites can thus be obtained by calculating the 
motion in the Galaxy backwards in time and looking for 
crossings of the orbit with the Galactic mid-plane. The 
times in the past at which the plane crossings occur yield 
kinematic estimates for the age of the DNS, while com- 
parison of the system's total center-of-mass velocity with 
the local Galactic rotational velocity at the birth sites 
yields the system's post-SN peculiar velocity. 

Secondly, the orbital semi-major axis and eccentric- 
ity right after the SN explosion forming pulsar B are 
determined by integrating the equations governing the 
evolution of the orbit due to gravitational wave emission 
backwards in time. The integration is performed for each 
Galactic plane crossing found from the Galactic motion 
calculations. Since each crossing yields a different kine- 
matic age, and thus a different endpoint of the reverse 
orbital evolution calculation, the post-SN parameters are 
a function of the considered Galactic plane crossing. 

Thirdly, the conservation laws of orbital energy and 
angular momentum are used to map the post-SN binary 
parameters to the pre-SN ones. The mapping depends on 
the kick velocity imparted to pulsar B at birth and re- 
sults in constraints on the magnitude and direction of the 
kick velocity, the mass of pulsar B's pre-SN helium star 
progenitor, and the pre-SN orbital separation. The pre- 
SN orbital eccentricity is assumed to be zero, as expected 
from strong tidal forces operating on the binary during 
the mass-transfer phase(s) responsible for spinning up 
pulsar A. The kick and pre- and post-SN binary parame- 
ters are then furthermore subjected to the requirements 
that they be consistent with the post-SN peculiar veloc- 
ity obtained from the Galactic motion calculations and 
with the observationally inferred spin-orbit misalignment 
angle of pulsar A. The latter constraint requires the ad- 
ditional assumption that tidal forces align the pre-SN ro- 
tational angular momentum vector of pulsar A with the 
pre-SN orbital angular momentum vector. 

Fourthly, the solutions of the conservation laws of or- 
bital energy and angular momentum are weighted ac- 
cording to the likelihood that they lead to the system's 
currently observed position and velocity in the Galaxy. 
A PDF of the admissible kick velocity and progenitor pa- 
rameters is then constructed by binning the solutions in 
a multi-dimensional "rectangular" grid. The maximum 
of the resulting PDF yields the most likely magnitude 
and direction of the kick velocity imparted to pulsar B 
at birth, mass of pulsar B's pre-SN helium star progen- 
itor, and pre-SN orbital separation. We also investigate 
the sensitivity of the PDF to the adopted assumptions by 
systematically varying the underlying assumptions, such 
as, e.g., the age and the magnitude of the transverse sys- 
temic velocity component. 



IV. KINEMATIC AGE AND HISTORY 
A. Galactic motion 

Tracing the motion of PSRJ073 7-3039 back in time 
requires the knowledge of both its present-day position 
and 3-dimensional space velocity. Unfortunately, no 
method is presently available to measure radial velocities 
of DNSs, so that the knowledge of their space velocity 
is limited to the component perpendicular to the line- 
of-sight (transverse velocity). At present, only an upper 
limit is available for the double pulsar: Vt < 30kms~^. 
The direction of the motion in the plane perpendicular 
to the linc-of-sight is still unknown. Therefore, as in Pa- 
per II, we explore the system's kinematic history in terms 
of two unknown parameters: (i) the radial component Vr 
of the present-day systemic velocity, and (ii) the orien- 
tation f2 of the transverse velocity in the plane perpen- 
dicular to the line-of-sight. Since the results presented in 
this paper do not depend on the exact definition of O, we 
refer the reader to Paper II for a more detailed discussion 
and definition. 

The motion of the system in the Galaxy is calculated 
with respect to a right-handed frame of reference OXYZ 
with origin at the Galactic center and with the XY-p\a,ne 
coinciding with the Galactic mid-plane. The X-axis is 
chosen such that the Sun is located in the XZ-pla,ne, 
and the positive direction of the K-axis such that it coin- 
cides with the direction of the Galactic rotational veloc- 
ity at the position of the Sun. The velocity components 
Vx, Vy, Vz with respect to this frame of reference are 
obtained from the radial and transverse velocity compo- 
nents Vr and Vt by a standard transformation of vector 
components, similar to the one adopted in § 2.2 of Pa- 
per II. We note that because PSR J0737-3039 is located 
only 20 pc below the Galactic mid-plane [53], the un- 
known radial velocity Vr mainly affects the determination 
of Vx and Vy ■ 

As in Paper II, we calculate the motion of the system 
in the Galaxy backwards in time using the Galactic po- 
tential of Carlberg & Innanen [28] with updated model 
parameters derived by Kuijken & Gilmore [29] [54]. Since 
the present-day kinematical constraints only provide an 
upper limit of 30 km on the transverse systemic veloc- 
ity and the calculation of the past orbit requires precise 
starting values for both the present-day position and the 
present-day velocity, we calculate the motion backwards 
in time for two specific velocities: Vt = lOkms^^ and 
Vt = 30kms~^. For each of these Vt, we integrate the 
equations of motion up to 100 Myr back in time for all 
possible U,.-values between -1500 kms~^ and 1500 kms~^ 
(in steps of lOkms"^), and i7- values between 0° and 360° 
(in steps of 5°). Although extremely high radial velocities 
in excess of ^ 1000 km s^^ may seem rather unlikely, they 
cannot be firmly excluded by the presently known obser- 
vational constraints (see Paper II). The possibility of a 
high space velocity is also supported by observations of 
high- velocity single radio pulsars such as PSRB2224-I-65 
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(with Vt > 800kms-i [30]) and PSRB1508+55 (with 
Vt = lOSSigg-'^kms"^ [31]). Moreover, as many of these 
high-velocity pulsars are expected to quickly escape the 
Galaxy, there is a strong bias against finding them. The 
fraction of radio pulsars formed with high space veloci- 
ties is therefore not well constrained in the current pulsar 
sample. 

In order to account for the possibility that PSR J0737- 
3039 has a large radial velocity and examine the effect 
of incorporating large radial velocities in the analysis, we 
weigh each considered radial velocity according to a pre- 
determined radial velocity distribution. For this purpose, 
we performed a population synthesis study of Galactic 
DNSs, including their kinematic evolution in the poten- 
tial of Carlbcrg & Innancn [28] . Theoretical radial veloc- 
ity distributions are then obtained by taking a snapshot 
of the population at the current epoch and determining 
the radial velocity for each DNS in the sample. The re- 
sulting PDFs are found to be represented well by Gaus- 
sian distributions with means of Okms~^ and velocity 
dispersions of 60-200 km s^^, depending on the adopted 
population synthesis assumptions (see § IV B for more de- 
tails). For comparison, wc also consider a uniform radial 
velocity distribution in which all radial velocities between 
-1500 km s~^ and 1500 km s~^ are equally probable [55]. 
It is clear that, contrary to what is stated by Piran & 
Shaviv [21], we do not assume the system to be moving 
almost exactly towards us with a very large radial veloc- 
ity Vr- For the other unknown parameter, the orientation 
angle O of the transverse velocity component in the plane 
perpendicular to the line-of-sight, we assume a uniform 
distribution between 0° and 360°. 

As noted in the previous section, it is reasonable to as- 
sume that PSR J0737-3039 was born close to the Galac- 
tic plane, so that the intersections of the past orbits with 
the Z = plane yield possible birth sites as functions 
of Vr and Q. The number of Galactic plane crossings 
associated with each Vr and f2, within the adopted age 
limit of 100 Myr, can be anywhere between 1 and 5. The 
times in the past at which the system crosses the Galactic 
plane furthermore yield kinematic estimates for the age 
of the DNS, while subtraction of the Galactic rotational 
velocity from the total systemic velocity at the birth sites 
yields the peculiar velocity right after the formation of 
pulsar B. 

Since the post-SN peculiar velocity represents one of 
the most important elements in constraining the kick ve- 
locity imparted to pulsar B at birth, it is worth taking 
a closer look at the range and distribution of the post- 
SN peculiar velocities derived from the Galactic motion 
calculations. In Fig. 1, the distribution of post-SN pe- 
culiar velocities is shown for a present-day transverse ve- 
locity component of 30kms~^, and kinematic age ranges 
of 0-100 Myr and 49-51 Myr. These probability distribu- 
tions are calculated considering weights according to (i) 
the probability that the system has a present-day radial 
velocity Vr (assumed to be distributed according to a 
uniform or a Gaussian distribution), (ii) the transverse 



velocity has a direction angle 17 (assumed to be uniformly 
distributed), and (iii) the system is found at its current 
position in the Galax;y (determined by the time the sys- 
tem spends near its current position divided by its kine- 
matic age, see § IV C for details). In particular, the prob- 
ability that, for a given pair of Vr and the post-SN 
peculiar velocity of the binary is equal to Vcm is given 

by 

P{VcM\Vr,n)0C J2 ^17-^^^(^'-'")' 

where N{Vr, O) is the number of Galactic plane crossings 
associated with Vr and f2, T{Vr, fi) is the time the system 
spends near its current position for the orbit associated 
with Vr and 17, and Tkin,i (K-, 17) is the kinematic age 
corresponding to the i-th plane crossing associated with 
Vr and O. The factor Ai {Vr,fl) is equal to 1 when the 
i-th plane crossing along the orbit associated with Vr and 
17 yields a post-SN peculiar velocity equal to Vcm, and 
equal to otherwise. The total probability that the post- 
SN peculiar velocity of the binary is equal to Vcm is then 
obtained by integrating Eq. (1) over all possible values 
of Vr and 17: 

PiVcM)0i [ [ P{VcM\Vr,n)P{Vr)P{n)dVrdn. 
J a JVr 

(2) 

Here P{Vr) is the probability distribution of the unknown 
radial velocity Vr^ and J-'(17) the probability distribution 
of the unknown proper motion direction 17 in the plane 
perpendicular to the line-of-sight. 

For ages of 0-100 Myr, post-SN peculiar velocities up 
to lOOkms"^ are likely for all four radial velocity dis- 
tributions. The highest post-SN peculiar velocities are 
found for the Gaussian radial velocity distributions with 
velocity dispersions of 130 km s^^ and 200kms~^ (Vcm 
values remain likely up to 200-300 km s~^) and the uni- 
form radial velocity distribution (Vcm values follow an 
almost flat distribution from 300kms~^ to 1500 km s~^). 
For ages of 49-51 Myr, the post-SN peculiar velocity dis- 
tributions all have a peak at 50kms~^. In the case of 
the uniform radial velocity distribution and the Gaus- 
sian distributions with velocity dispersions of 130kms~^ 
or 200kms~^, additional peaks of almost equal height 
occur at even higher post-SN peculiar velocities. Similar 
conclusions apply to the post-SN peculiar velocity dis- 
tributions obtained for a present-day transverse velocity 
component of lOkms"^. Despite the small lower limit 
on the present-day transverse systemic velocity, a wide 
range of non-negligible post-SN peculiar velocities is thus 
possible. We note however that the distributions incor- 
porate all possible post-SN peculiar velocities obtained 
from tracing the motion of the system backwards in time 
and that some of these may require SN kicks and mass 
loss that are incompatible with the observational con- 
straints on the orbital elements and component masses 
of the double pulsar. 
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FIG. 1: Distribution of post-SN peculiar velocities for a present-day transverse velocity of SOkrns^^, kinematic ages ranges 
of 0-100 Myr and 49-51 Myr, and radial velocity distributions varying from a uniform distribution (solid line) to Gaussian 
distributions with velocity dispersions of 60kms~^ (long-dashed line), 130kms~^ (short-dashed line), and 200kms~^ (dash- 
dotted line). For clarity, the PDFs are offset from each other by an arbitrary amount. 
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FIG. 2: Distribution of kinematic ages for a present-day transverse velocity of SOkms^"'^, kinematic age ranges of 100 Myr 
and 30-70 Myr, and radial velocity distributions varying from a uniform distribution (solid line) to Gaussian distributions with 
velocity dispersions of 60kms~^ (long-dashed line), 130kms~^ (short-dashed line), and 200kms~^ (dash-dotted line). For 
clarity, the PDFs are ofTset from each other by an arbitrary amount. 



Similarly, PDFs can be constructed for the kinematic 
ages Tkin- They are shown in Fig. 2 for a present-day 
transverse velocity component of 30kms~^, and kine- 
matic age ranges of 0-100 Myr and 30-70 Myr. For ages 
of 0-100 Myr, the PDFs are all very strongly peaked at 
very young ages. In particular, 90-97% of the Galactic 
plane crossings give rise to ages of 0-5 Myr. The peak of 
the PDF at ages of 0-5 Myr furthermore increases with 
increasing velocity dispersion of the adopted radial veloc- 
ity distribution. For ages of 30-70 Myr, the PDFs over- 
all decrease with increasing values of ridn and show no 
strong preference for any particular range of kinematic 
age values. Similar conclusions apply to the kinematic 
ages obtained for a present-day transverse velocity com- 
ponent of lOkms"^. 



B. Theoretical radial velocity distributions 



While the radial velocity of PSR J0737-3039 is un- 
known, population synthesis calculations can provide us 
with some reasonable guidance on the expected radial 
velocity distribution of the DNS population. We use the 
StarTrack population synthesis code (see Belczynski et 
al. [9]; and especially Belczynski et al. [32] for an ex- 
tensive and detailed description of this binary popula- 
tion synthesis evolutionary code) to develop models of 
DNS formation and evolution in our Galaxy (assuming 
continuous star- formation history over lOGyr and solar 
metallicity). To calculate the DNS kinematic evolution 
we need to also choose a spatial distribution of DNS bi- 
nary progenitors in the Galaxy. We adopt a standard 
double exponential to represent the density distribution 
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FIG. 3: Present-day radial-velocity distributions of coalescing DNSs for a population synthesis model with Maxwellian kicks 
with a velocity dispersion Ukick = 250kms~^. The left panel shows the distribution for the entire population of coalescing 
DNSs, while the right panel shows the distribution of coalescing DNSs within 50 pc of the Galactic mid-plane and with a proper 
motion of less than 30kms~^. The dashed lines in the left and right panels are best-fit Gaussian distributions with velocity 
dispersions of 117kms~^ and ISOkms"^, respectively. Note that in the right panel 45% of the systems have \Vr\ > 100 kms"'^. 



TABLE II: Velocity Dispersions for Supernova Kicks and 
Radial Velocities 



CTkick 


CVr'' 


(kms-1) 


(kms-1) 


50 


63 ±2 


100 


94 ±4 


150 


92 ±4 


200 


104 ±5 


250 


130 ±7 


300 


148 ±9 


350 


124 ±7 


400 


172 ± 15 


450 


147 ± 14 


500 


280 ± 141 



"'Maxwellian distribution. 

'Best fit Gaussian distribution centered at zero. Values are for 
coalescing DNSs within 50 pc of the Galactic mid- plane and with 
a proper motion of less than 30kms~^. 

of Galactic systems: 

n{R, z) = no exp {-R/Hr) exp {-\z\/Hz) , (3) 

where rio = 1/ {AttHzH'^^ is the normalization, Hn the 
scale length, and Hz the scale height. We adopt Hr = 
2.8 kpc and Hz — 70 pc, appropriate for massive stars. 
The motion of the DNSs in the Galaxy is calculated using 
the same Galactic potential as in §IVA. 

In this analysis, we consider a set of models with 
input parameters chosen as for the standard model 
(model A) in Belczynski et al. [9], except for the as- 
sumed Maxwellian distributions of the kick magnitudes 
imparted to nascent neutron stars: we vary the one- 
dimensional kick velocity dispersion (Ttick in the range 
of 50 — 500 km s"-*^, as shown in Table II. We restrict 
our parameter study to these models, because a broad 



model exploration is not the primary goal of our study 
and because prior experience with DNS models clearly 
indicates that the assumed NS kick velocity model is the 
most important parameter affecting the DNS kinematic 
evolution and properties [12]. 

We use this set of DNS models to derive radial-velocity 

distributions with respect to the Sun. To obtain a realis- 
tic representation of the kinematic properties of the cur- 
rent population of tight DNS systems (like PSRJ0737- 
3039), we place a cut on the model DNS populations and 
select only those that have already formed (i.e., they have 
a formation time shorter than lOGyr) and will coalesce 
within the age of the Universe (at a time in the range 
9-14 Gyr). We next distribute the selected DNS popu- 
lation in the Galaxy and follow their motion using the 
calculated center-of-mass velocities imparted to the bi- 
naries after both of the two supernova events [56]. We 
also take into account the effect of Galactic rotation. Wc 
then derive the radial velocity distribution by extract- 
ing the calculated peculiar velocities (relative to the Sun) 
and considering the appropriate projections. We can fur- 
ther impose a number of additional constraints on the 
model population to select systems that have kinematic 
properties similar to PSR J0737-3039: coalescing DNS 
with proper motion smaller than 30kms~^ and a loca- 
tion within 50 pc from the Galactic plane; we note that 
we do not impose any constraints on the binary proper- 
ties to avoid any statistical accuracy problems. 

We find that the derived radial velocity distributions 
are very well approximated by a Gaussian centered 
at zero with velocity dispersions ctv, ranging from 60 
to 200 km s"-'^. The dependence of the derived radial- 
velocity dispersions on the kick velocity dispersions are 
shown in Table II. The distributions with and with- 
out the kinematic constraints arc shown in Figure 3, for 
Ckick = 250kms~^. The corresponding values of ctVp 
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for the two panels are rather similar: 117kms^^ and 
130kms~^, respectively. Moreover, for this akick value, 
45% of the systems have \Vr \ > lOOkms"^. 

It is important to note that imposing the proper- 
motion constraint does not actually reduce the typical 
magnitude of the radial velocity. The reason is that small 
proper motion values select two sub-populations of DNS 
systems: one with inherently small center-of-mass veloc- 
ities (and consequently small radial velocities), and an- 
other with large center-of-mass velocities but appropri- 
ate direction so that the proper motion is small. This 
second group naturally leads to large radial velocities 
because, for a given center-of-mass velocity, the projec- 
tion affects the radial component in the opposite way 
as it does the transverse component. Piran & Shaviv 
[21] somewhat mislcadingly draw the focus away from 
this second group and use the first sub-population to ar- 
gue that both the present and post-SN peculiar velocity 
must be small. Since our theoretical radial-velocity dis- 
tributions incorporate the constraints \Z\ < 50 pc and 
Vt < 30kms~^, they implicitely account for the geomet- 
rical probability argument by Piran & Shaviv [21] that 
the probability of a small transverse velocity component 
is larger for smaller total systemic peculiar velocities (see 
their Eq. 2). Moreover, as outlined in detail in the pre- 
vious section, we explicitely account for differences be- 
tween the system's current and post-SN peculiar velocity 
by modeling its past motion in the Galactic potential. 

Last, we note (although not directly relevant to this 
study) that with the full population synthesis models 
we confirm the results of Kalogera et al. [33] about the 
vertical DNS distribution: for (Jkick = 250kms~^, it is 
best fitted by two exponentials, each with a typical scale 
height of 0.230 ± 0.006 kpc and 1.72 ± 0.07 kpc. As ex- 
pected the scale height increases with increasing typical 
kick magnitude. 



C. Probability to find the system at its current 
position 

At a distance of ~ 600 pc and a Galactic latitude of 
-4.5°, PSR J0737-3039 is currently located only 20 pc 
below the Galactic plane. Piran & Shaviv [19] used the 
proximity of the system to the Galactic plane to argue 
that pulsar B most likely received only a small kick ve- 
locity at birth and that its helium star progenitor most 
likely had a pre-SN mass only slightly larger than the 
pulsar's present mass. Their procedure, however, only 
accounted for the motion of the system perpendicular to 
the Galactic plane. 

Here, we determine the probability of finding the sys- 
tem at its current position in the Galaxy in three dimen- 
sions (i.e., besides the vertical distance to the Galactic 
plane, we also consider the radial and azimuthal position 
in the plane). For this purpose, we determine the time 
the system spends in a sphere with radius R centered on 
its current location for all Vr- and 0-values considered 



for the derivation of the possible birth sites in §IVA. 
The probability to find the system near its current po- 
sition is then determined as the time it spends in this 
sphere divided by its age. For the latter, we adopt the 
kinematic ages obtained by tracing the motion of the sys- 
tem backwards in time, so that the probability of finding 
the system near its ciurent position depends on the ra- 
dial velocity Vr, the proper motion direction fl, and the 
Galactic plane crossing considered along the orbit associ- 
ated with Vr and f2. We performed some test calculations 
to successfully verify that the results based on the thus 
determined probabilities are insensitive to the adopted 
value of R, as long as it is sufficiently small for the sphere 
to represent a local neighborhood near the system's cur- 
rent position. For the results presented in this paper we, 
somewhat arbitrarily, choose R = 50 pc. 



V. PROGENITOR CONSTRAINTS 

A. Post-supernova orbital peirameters 

After the formation of pulsar B, the evolution of the 
system is expected to be driven exclusively by the emis- 
sion of gravitational waves. The post-SN orbital semi- 
major axis and eccentricity can then obtained from the 
current values by numerically integrating the system of 
differential equations governing the orbital evolution due 
to gravitational wave emission backwards in time. For 
this purpose, we use the differential equations derived by 
Junker & Schafer [34] which are valid up to 3.5 post- 
Newtonian order of approximation. 

We calculate the orbital evolution backwards in time 
for all Galactic plane crossings associated with all values 
of Vr and Q. considered in the calculation of the possible 
birth sites. For each crossing the integration is termi- 
nated at the associated kinematic age Tkin- We find the 
resulting post-SN orbital separation A to be always be- 
tween 1.26 i?0 and 1.54 i?©, and the post-SN orbital ec- 
centricity e between 0.088 and 0.12 (see Papers I and II). 
Since the kinematic ages are functions of the radial veloc- 
ity Vr, the proper motion direction f2, and the Galactic 
plane crossing considered along the past orbit associated 
with Vr and fi, the particular values of A and e associated 
with a given Tkin are also functions of these variables. 



B. Orbital dynamics of eisymmetric supernova 
explosions 

The pre- and post-SN binary parameters and the kick 
velocity imparted to pulsar B at birth are related by the 

conservation laws of orbital energy and angular momen- 
tum. For a circular pre-SN orbit, the relations take the 
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form [e.g. 35-37, 37-39] 

V^ + V^ + 2V^Vo cose 

= G(M. + M.)(-1_^ 

Al [V^ sin^ e cos^ (p + (14 cos 61 + ^0)^ 
= G{MA + MB)A{l-e''), 



(4) 



(5) 



where Mq is the mass of pulsar B's pre-SN helium 

star progenitor, Aq the pre-SN orbital separation, Vq = 
[G{Ma + Mo)/^o]^/^ the relative orbital velocity of pul- 
sar B's pre-SN helium star progenitor, and 14 the mag- 
nitude of the kick velocity imparted to pulsar B at birth. 
The angles 9 and 4> define the direction of the kick ve- 
locity imparted to pulsar B: 6 is the polar angle between 
Vk and Vo, and (f> the corresponding azimuthal angle in 
the plane perpendicular to Vq defined so that (/) = 7r/2 
corresponds to the direction from pulsar A to pulsar B's 
helium star progenitor (see Fig. 1 in [22] for a graphical 
representation) . 

The requirement that Eqs. (4) and (5) permit real so- 
lutions for Vk, 0, (j), Mo, and Aq, imposes constraints 
on the pre- and post-SN binary parameters and on the 
magnitude and direction of the kick velocity imparted to 
pulsar B at birth. For a mathematical formulation of 
these constraints, we refer to Eqs. (21)-(27) in [40], and 
references therein (a more compact description can also 
be found in Papers I and II). We here merely recall that 
the constraints express that: (i) the binary components 
must remain bound after the SN explosion, (ii) the pre- 
and post-SN orbits must pass through the instantaneous 
position of the component stars at the time of the SN 
explosion, and (iii) there is a lower and upper limit on 
the degree of orbital contraction or expansion associated 
with a given amount of mass loss and a given SN kick. 

Besides the changes in the orbital parameters and the 
mass of the exploding star, the SN explosion also imparts 
a kick velocity to the binary's center of mass and tilts the 
post-SN orbital plane with respect to the pre-SN one. 
Under the assumption that the pre-SN peculiar velocity 
of the binary is small in comparison to the pre-SN orbital 
velocity of the component stars, the magnitude of the 
post-SN center-of-mass velocity Vcm is given by [e.g. 37] 



•'cm — 



Mb Mo 



+ 



{Ma + Ma){MA + Mb) 



G{Mo - Mb)Ma 



MqA 



G{Mo - Mb){Mo - 2 Mb)Ma 
Mo Mb Ao 



(6) 



The tilt angle A between the pre- and post-SN orbital 
planes is given by [22] 



cos A 



A Ma -\- Mb 
To Ma + Mo 



Vo 



(1-e^) 



-1/2 



(7) 



From Eqs. (5) and (7), one derives that A = 0° for kicks 
directed in the pre-SN orbital plane (i.e. cos9 — ±1), 
independent of the kick velocity magnitude 14 . 



C. Stellar and binary evolution 

The constraints on the progenitor of PSR J0737-3039 
resulting from the dynamics of asymmetric SN explosions 
arise solely from tracing the evolution of the current sys- 
tem properties backwards in time. The pre-SN orbital 
separations and pulsar B progenitor masses found this 
way are, however, not necessarily accessible through the 
currently known DNS formation channels. Further con- 
straints on the progenitor of pulsar B can therefore be 
obtained from stellar and binary evolution calculations. 

A lower limit on the mass of pulsar B's pre-SN pro- 
genitor is given by the requirement that the star must 
be massive enough to evolve into a NS rather than a 
white dwarf. According to our current understanding of 
helium star evolution, the minimum helium star mass 
required for NS formation is 2.1-2.8Mq [10, 41]. The 
actual helium star minimum mass is, however, still con- 
siderably uncertain due to the poorly understood evolu- 
tion of massive stars and possible interactions with close 
binary companions. 

As in Papers I and II, we impose a lower limit of 2.1 Mq 
on the mass of pulsar B's pre-SN helium star progeni- 
tor. However, in the light of recent suggestions that the 
progenitor of pulsar B may have been significantly less 
massive than the conventional lower limit of 2.1 Mq, we 
also explore the possibility of progenitor masses as low 
as pulsar B's present-day mass of 1.25 [3, 19]. Un- 
less our current understanding of helium star evolution 
is seriously flawed, this scenario implies that the progen- 
itor of pulsar B must have lost a significant amount of 
mass (at least ~ 0.7 M©) after it had already established 
a sufficiently massive core to guarantee the occurrence 
of a future SN explosion. We note that this possibility 
is included in the binary population synthesis calcula- 
tions used to derive the adopted theoretical DNS radial- 
velocity distributions (see Table II). Contrary to what is 
stated by Piran & Shaviv [21], our population synthesis 
calculations do allow for low-mass helium-rich progen- 
itors of the second NS and the derived radial-velocity 
distributions are therefore by no means a priori biased 
towards high progenitor masses. 

In Paper I, we also showed that the pre-SN binary 
was so tight (1.2 i?© < < 1.1 Rq) that the helium 
star progenitor of pulsar B must have been overflowing 
its Roche lobe at the time of its SN explosion (see also 
[13]). An upper limit on the progenitor mass is therefore 
given by the requirement that this mass-transfer phase 
be dynamically stable (otherwise the components would 
have merged and no DNS would have formed) [57]. Based 
on the evolutionary tracks for NS -|- helium star bina- 
ries calculated by Ivanova et al. [42], we adopt an up- 
per limit of 3.5 on the mass ratio Mq/Ma of the pre- 
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SN binary to separate dynamically stable from dynam- 
ically unstable Roche- lobe overflow (see also [43, 44]). 
For Ma = 1.34 M©, this yields an upper limit of 4.7 
on the mass of pulsar B's pre-SN helium star progenitor. 
Because of the mass transfer, the pre-SN progenitor mass 
may, however, be considerably lower than the progenitor 
mass at the onset of Roche-lobe overflow. The upper 
limit of 4.7M0 is therefore fairly conservative. However, 
as we will see in the next section, the preferred pre-SN 
progenitor masses are always significantly smaller than 
4.7 Mq so that a stricter upper limit would not consid- 
erably affect any of our results. 



VI. PROBABILITY DISTRIBUTION 
FUNCTIONS 

A. A coherent picture for the progenitor and 
formation of pulsar B 

The steps outlined in the previous sections have been 
used to derive the most likely pre-SN progenitor of 
PSR J0737-3039 and the magnitude and direction of the 
kick velocity imparted to pulsar B at birth. Firstly, the 
possible kinematic ages and post-SN peculiar velocities 
are determined by tracing the Galactic motion backwards 
in time using the current constraints on the proper mo- 
tion and position in the Galaxy. These ages and veloci- 
ties are functions of the presently unknown values of the 
radial velocity Vr and proper motion direction il. Sec- 
ondly, the kinematic ages are used to reverse the orbital 
evolution due to gravitational radiation and determine 
the orbital semi-major axis and eccentricity right after 
the formation of pulsar B. Thirdly, the pre-SN binary pa- 
rameters and kick velocity imparted to pulsar B are con- 
strained using the conservation of orbital energy and an- 
gular momentum during asymmetric SN explosions. The 
constraints are supplemented with the requirements that 
the post-SN peculiar velocity obtained from Eq. (6) be 
compatible with the post-SN peculiar velocity obtained 
by tracing the motion of the system in the Galaxy back- 
wards in time, and that the spin-orbit misalignment of 
pulsar A obtained from Eq. (7) be compatible with the 
observational constraint |A — 90° | > 30°. Imposing the 
latter constraint rests on the additional assumption that 
pulsar A's spin is aligned with the pre-SN orbital angu- 
lar momentum, as expected from the strong tidal inter- 
actions operating during the spin-up phase of pulsar A. 



B. The most likely kick velocity and progenitor 
properties 

Since the conservation laws of orbital energy and an- 
gular momentum yield two equations for five unknown 
quantities Vfe, cos^, (j), Mq, and Aq, we solve the equa- 
tions for Mo and as functions of Vfe, cos 6', and 0. 
Each solution is then weighted assuming uniform prior 



TABLE 111: Effect of marginalizing the 5-D PDF for Vfc, cos 6», 
(ji, Mo, and Aq on the most likely kick velocity and progenitor 
parameters for a transverse velocity component of 30kms~^, 
a minimum helium star mass of 1.25 Mq, an age range of 
0- 100 Myr. and a uniforiri prcsciit-day radial velocity distri- 
bution. 







Most likely value 




Variables 


Vk 


cos 6 


4> 


Mo 


Ao 




(kms"'; 




(deg) 


(Mo) 


(Re) 


Vk, cos 61, 0, Mo, Ao 


75 


-0.025 


10 


1.45 


1.15 


Vk, cos 9, 0, Mo 


75 


-0.025 


10 


1.45 




Vk, cose, Mo 


45 


-0.050 




1.45 




Vk, Mo 


10 






1.45 




Vk 


55 











distributions for Vk, cos 6, (f), Mq, and ^o- Thus the 
kick direction is assumed to be isotropically distributed 
in space. Similar to Eqs. (1) and (2), we next determine 
the probability that, given the currently known obser- 
vational constraints, PSR J0737-3039B was born with a 
natal kick velocity Vk with a direction determined by 
cos 9 and (f> from a progenitor of mass Mq orbiting the 
first-born pulsar in a circular orbit of radius Aq as 



P{Vk,cos0, (l),Mo,Ao\Xobs) 



oc 



n JVr- 



J2 3l%2L„.,v,.n) 



xP{Vr)P{n)dVrdn. 



(8) 



Here ATobs denotes all currently known observational con- 
straints (orbital separation, eccentricity, proper motion, 
etc). The factor Kj {Vr,fl) is equal to 1 when the con- 
sidered Vfe, cos 6', (p, Mq, and Aq satisfy all constraints 
imposed by tracing the system's kinematic and evolu- 
tionary history back to the formation time of pulsar B 
for the i-th plane crossing along the orbit in the Galaxy 
associated with K and fl, and equal to otherwise. The 
other quantities are defined in § IV A. Hence, we obtain 
a 5-dimensional PDF for the kick velocity and progenitor 
parameters: Vk, cos 9, cj), Mq, and Aq. 

The construction of the full 5-D PDF represents a 
significant departure from previous investigations which 
restricted the analysis to the derivation of 1-D PDFs 
for Vfe (Papers I and II) or 2-D PDFs for Vfe and Mq 
[19, 20]. Although the derivation of the most likely Vfe 
and Mq values from these marginalized PDFs is math- 
ematically correct, they may be subjected to projec- 
tion effects which can drastically affect the position of 
the maxima. In order to illustrate this, we constructed 
the 5-D PDF for a transverse velocity component of 
30kms~^, a minimum helium star mass of 1.25 M0, an 
age range of 100 Myr, and a uniform present-day ra- 
dial velocity distribution. The most likely kick velocity 
and progenitor properties obtained from the 5-D PDF are 
Vfe ~ 75kms-i, cos 61 ~ -0.025, = 10°, Mq ~ 1.45 Mq, 
and Aq ~ 1.15 i?0. The most likely values obtained af- 
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ter successively integrating over Aq, (j), cos 9, and Mq are 
listed in Table III. While the integration over does not 
significantly affect the most likely value of Vk , successive 
integrations over <p and cos 9 reduce the most likely Vk to 
45kms~^ and lOkms"^, respectively. A subsequent in- 
tegration over Mo furthermore increases the most likely 
Vk back to 55kms~^. Marginahzing a multi-dimensional 
PDF can thus have a significant impact on the location 
of its peaks. Some caution is therefore in order when 
interpreting marginalized PDFs. 

The most likely kick velocity and progenitor param- 
eters can also depend strongly on the assumptions un- 
derlying the reconstruction of the binary's evolutionary 
history (i.e., the magnitude of the transverse systemic 
velocity component, the minimum helium star mass re- 
quired for the formation of a NS, the age of the system, 
and the probability distribution of the unknown radial 
velocity). In order to investigate this, we systematically 
calculated the 5-D PDFs for Vk, cos9, (f), Mq, and Aq for 
different sets of prior assumptions. The resulting most 
likely kick velocity and progenitor parameters are sum- 
marized in Table IV. The dependence of the most likely 
Vk, con 9, and Afo on the adopted prior assumptions is 
also presented graphically in Figs. 4-6. For brevity, we 
denote the uniform radial velocity distribution in the ta- 
ble and figures as a Gaussian distribution with an infinite 
radial velocity dispersion. 

The most likely kick velocity imparted to pulsar B at 

birth is smaller than 50kms~^ only when the transverse 
velocity component has a magnitude of lOkms"^, the 
minimum helium star mass required for NS formation is 
allowed to be as low as 1.25 Mq, and the age of the system 
is assumed to be between 30 and 70 Myr. All other models 
yield most likely kick velocities of 50-180 km s^ ^ . When a 
minimum pre-SN helium star mass of 2.1 Mq is imposed, 
the kicks are always strongly favored to be directed op- 
posite to the helium star's pre-SN orbital motion (most 
likely cos9 ~ — 0.90 ± 0.05). When the constraint on the 
minimum pre-SN helium star mass is relaxed, the most 
likely kick direction can shift significantly and can even 
become perpendicular to the helium star's pre-SN orbital 
velocity. Allowing pre-SN helium star masses down to 
1.25^^0 furthermore always leads to most likely progen- 
itor masses of 1.3 — 2.0 Mq, except when Vt = lOkms"^, 
the age of the system is between and 100 Myr, and the 
radial velocities are distributed uniformly or according to 
a Gaussian with a velocity dispersion of 200kms~^. 

In summary, while small kick velocities of just a few 
tens of kms~^ could be favored for some models, the 
majority of the models yields most likely values of 50- 
180 kms"-"^. Progenitor masses below 2.1 Mq are further- 
more not required to explain the system properties, al- 
though they arc generally favored when helium stars be- 
low 2.1 Mq are still assumed to be viable NS progenitors. 
We note though that the results presented in Table IV 
and the associated figures are based on the assumption 
that all kick velocity magnitudes Vk are equally proba- 
ble. This assumption is inconsistent for models using the 



Gaussian radial velocity distributions based on popula- 
tion synthesis calculations of coalescing DNSs. In par- 
ticular, these calculations all adopt a MaxwcUian rather 
than a uniform kick velocity distribution. In all cases, 
weighing the kick velocities according to the Maxwellian 
underlying the derivation of the radial velocity distribu- 
tions would, however, shift the most likely kick velocities 
to higher Vk values than listed in Table IV. This rein- 
forces our conclusion that the presently known observa- 
tional constraints not necessarily disfavor kick velocity 
magnitudes of lOOkms"^ or more. 



C. Confidence limits on the kick velocity and 
progenitor pEirameters 

The use of the full 5-D PDF for Vk, cos 6*, (j), Mq, and 
^0 in the derivation of the most likely pulsar B kick ve- 
locity and progenitor parameters has the advantage of 
being free of projection effects, but is much more cum- 
bersome and computationally expensive for the calcula- 
tion of confidence limits. Since our main aim in this 
paper is to present a differential analysis showing the 
effects of different model assumptions, we adopt a sim- 
pler approach and determine the confidence limits from 
marginalized 1-D PDFs. We recognize, however, that 
in doing so, projection effects may play a significant role 
and any information on possible correlations between the 
different parameters will be lost. 

For illustration, some representative 1-D PDFs used 
for the calculation of the confidence limits in the case of 
a present-day transverse velocity component of 10 km s^^ 
are shown in Fig. 7. Panels (a)-(c) show the kick velocity 
distributions resulting from different present-day radial 
velocity distributions for ages ranges of 0-100 Myr and 
49-51 Myr, and minimum pre-SN helium star masses of 
1.25 Mq and 2.1 M0. For a given age range and mini- 
mum pre-SN helium star mass, the PDFs show a peak 
which is most pronounced when a Gaussian radial veloc- 
ity distribution with a velocity dispersion of 60kms^^ 
is considered, and which widens with increasing radial 
velocity dispersions. When the age is assumed to be 
0-100 Myr and the minimum pre-SN helium star mass 
1.25 Mq, there is a clear tendency for the 1-D PDFs to 
favor kick velocities of 50 km s~^ or less (although this be- 
comes significantly less pronounced with increasing radial 
velocity dispersions). This trend shifts towards favoring 
kick velocities of 50-100 km s~^ when the age range is 
narrowed to 49-51 Myr. Moreover, when the age range 
is kept fixed at 0-100 Myr, but the minimum pre-SN he- 
lium star mass is increased to 2.1 Mq, the favored range 
of kick velocities shifts to 100-150 km s~^. In the latter 
case, the kick velocity is furthermore always required to 
be larger than ^ GOkms"-'^ (see also Papers I and II). 

In Panels (d)-(e), we show the distribution of cos 9 for 
different present-day radial velocity distributions, DNS 
ages of 0-100 Myr, and minimum pre-SN helium star 
masses of 1.25 Mq and 2.1 Mq. For a minimum he- 
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TABLE IV: The most likely kick velocity and progenitor parameters obtained from the 5-D PDF for Vk, cosO, 4>, Mo, and 
Aq, for different transverse systemic velocity components, minimum pre-SN helium star masses, DNS ages, and present-day 
radial velocity distributions. For brevity, the uniform radial velocity distribution is denoted by av^ = oo. Tabulated values 

of (j) are restricted to 0° < < 90°. For each listed value additional and equally likely solutions to Eqs. (4) and (5) are 
associated with — </> and 180° ±0. For 14 and Mo, 68% and 90% confidence intervals are also listed. The intervals are obtained 
by marginalizing the 5-D PDF to 1-D PDFs for Vk and Mq. 



Priors Most likely progenitor 68% confidence intervals 90% confidence intervals 
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Hum star mass of 1.25 and Gaussian radial veloc- 
ity distributions with velocity dispersions of more than 
130kms~^, the PDFs favor kick directions opposite to 
the pre-SN orbital motion. In the case of the Gaussian 
radial velocity distribution with a velocity dispersion of 
60kms~^, the PDF behaves qualitatively entirely differ- 
ent and shows a slight preference for kicks perpendicu- 
lar to the pre-SN orbital velocity. A non-negligible low- 
probability tail furthermore extends well into the region 
of kick directions in the sense of the pre-SN orbital mo- 
tion. When the minimum mass for pulsar B's pre-SN 
helium star progenitor is raised to 2.1 Mq, the kicks are 
always directed opposite to the orbital motion (a lower 



limit on of ~ 113° was already derived in Papers I 

and 11) . The PDF furthermore becomes fairly insensitive 
to the adopted radial velocity distribution. 

Panel (f), finally, shows the distribution of possible 
pre-SN progenitor masses for different present-day ra- 
dial velocity distributions. DNS ages of 0-100 Myr, and 
a minimum pre-SN helium star mass of 1.25 M©. The 
distributions all favor progenitor masses of 1.4-1.5 M©, 
with the preference for this mass range being strongest 
for present-day radial velocity distributions with small 
radial velocity dispersions. Distribution functions for a 
minimum pre-SN helium star mass of 2.1 Mq look prac- 
tically the same as the ones displayed panel (f) if they 
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FIG. 4: Most likely SN kick velocity magnitude 14 for different sets of prior assumptions, for three age ranges: 0-100 Myr (solid 
lines), 30- 70 Myr (short-dashed lines), 49-51 Myr (long-dashed lines). For brevity, the uniform radial velocity distribution is 
denoted by avv = oo. 



were cut off at 2.1 Mq. 

Adopting a present-day transverse velocity of 
30kms^^ instead of lOkms"'^ yields very similar conclu- 
sions. Overall the effects of varying model assumptions 
tend to be somewhat less pronounced though. 

The 68% and 90% confidence limits on Vk and Mq ob- 
tained from the 1-D PDFs are listed in Table IV next to 
the most likely kick velocity and progenitor parameters 
[58]. In the majority of the cases, the confidence lim- 
its contain the most likely Vk and Mq obtained from the 
full 5-D PDF for Vk, cos 6*, 0, Mq, and Aq. Whenever 
this is not the case, projections effects cause significant 
shifts between the values favored by the 1-D and 5-D 
PDFs. The confidence limits on Vk readily show that the 
possibility of small kick velocities is entirely due to the 
assumption that helium stars less massive than 2.1 Mq 
can be viable NS progenitors. Under this assumption, 
the upper limits on 14 at the 68% confidence level vary 
from 65kms~^ to 275kms~^. Thus, even though small 
kick velocities become possible, larger kick velocities in 
excess of lOOkms"^ are all but statistically disfavored. 
This is even more true when the upper limits at the 90% 
confidence level arc considered. For a minimum prc-SN 
helium star mass of 2.1 Mq, the confidence limits for Vk 
shift to significantly higher values which are conform with 
the kick velocity magnitudes derived from observations of 
single radio pulsar [46-48]. These "conventional" kicks 



are also consistent with the tilt between pulsar B's rota- 
tion axis and the normal to the orbital plane inferred by 
Lyutikov & Thompson [49]. When the minimum pre-SN 
helium star mass is assumed to be 1.25 Mq, the 68% 
and 90% confidence limits on Mq typically go all the 
way down to this limit, but often also extend upwards to 
masses larger than 2.1 Mq (especially at the 90% confi- 
dence level) . These numbers confirm our previous conclu- 
sion that the presently known observational constraints 
not necessarily disfavor progenitors more massive than 
2.1 Mq and kick velocities in excess of lOOkms"^. 



D. Birth sites 

Just like the kick velocity and progenitor parame- 
ters, the most likely birth site of PSR J0737-3039B de- 
pends strongly on the adopted transverse systemic veloc- 
ity component, age range, minimum helium star mass, 
and present-day radial velocity distribution. In § IV A, 
we have shown that the distribution of kinematic ages 
in the range from to 100 Myr peaks strongly at ages 
of 1-2 Myr (sec Fig. 2). This implies that the system 
can travel only a short distance (less than 2kpc if the 
systemic velocity is less than 1000 km s~^) before reach- 
ing its current position. Adopting an age range of - 
100 Myr therefore yields PDFs for the birth sites that 
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FIG. 5: Most likely SN kick direction cos 6 for different sets of prior assumptions, for three age ranges: 0-100 Myr (solid lines), 
30-70 Myr (short-dashed lines), 49-51 Myr (long-dashed lines). For brevity, the uniform radial velocity distribution is denoted 
by = oo. 



are strongly peaked near the system's current position at 
X = — 0.5kpc and Y = — 8.3kpc, regardless of any of the 
other model assumptions. The situation changes dras- 
tically for the age ranges of 30-70 Myr and 49-51 Myr, 
which exclude very young ages and therefore impose a 
minimum distance traveled by the system since the time 
it was born. 



For illustration, the distribution of the possible birth 
sites in the Galactic plane is shown in Fig. 8 in the 
case of a prcscnt-day transverse velocity component 
of 30kms~^, a minimum pre-SN helium star mass of 
1.25 Mq, an age range of 30-70 Myr, and a Gaussian 
radial velocity distribution with a velocity dispersion of 
130 km s^^. The most likely birth sites are all located 
within a radius of 2kpc from {X,Y) = (—3.5, — 7.5)kpc. 
None of the birth sites are furthermore close to the sys- 
tem's present-day position indicated by the solid circle at 
{X,Y) = (—0.5, — 8.3)kpc. If the minimum helium star 
mass required for NS formation is assumed to be 2.1 Mq 
instead of 1.25 Mq, the most likely birth sites shown in 
Fig. 8 become entirely inaccessible. Instead, the binary 
then most likely originates from a circle with a radius of 
1 kpc centered on {X, Y) = (—2.5, —3.5) kpc, with a clear 
non-negligible lower probability tail extending more or 
less linearly towards {X,Y) = (3.5, 0.5) kpc. 



VII. COMPARISON WITH PREVIOUS WORK 

A. Willems et al. [15] 

In Paper II, we derived constraints on the formation 
and progenitor of PSR J0737-3039B, using the then avail- 
able observational constraints on the binary's properties 
and motion in the Galaxy. The main differences with the 
observational and theoretical information adopted in this 
investigation are the lower transverse systemic velocity, 
the constraint on pulsar A's spin-orbit misalignment an- 
gle, the possibility of a significantly lower minimum mass 
for pulsar B's pre-SN helium star progenitor, and the use 
of DNS population synthesis calculations to obtain the- 
oretical present-day radial-velocity distributions. More- 
over, in Paper II we mainly focused on the ranges of pos- 
sible progenitor parameters, and restricted the statistical 
analysis to the derivation of 1-D PDFs for pulsar B's na- 
tal kick velocity and pulsar A's spin-orbit misalignment. 
Among the observational and theoretical improvements, 
the availability of the constraint on pulsar A's spin-orbit 
misalignment has the least effect on the derivation of the 
possible pulsar B progenitor and kick velocity proper- 
ties. The main reason for this is that the excluded range 
of tilt angles from 60° to 120° was already statistically 
disfavored by the kick velocity and progenitor constraints 
derived in Paper II (see Fig. 12 in that paper). 
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FIG. 6: Most likely pro-SN progenitor mass Mo for different sets of prior assumptions, for three age ranges: 0-100 Myr (solid 
lines), 30-70 Myr (short-dashed lines), 49-51 Myr (long-dashed lines). For brevity, the uniform radial velocity distribution is 
denoted by av^ = oo. 



While the majority of the results presented in Paper II 
were obtained assuming a minimum pre-SN helium star 
mass of 2.1 Mq, we also discussed some preliminary test 
calculations allowing pre-SN helium star masses below 
2.1 Mq. Because of the strong influence on the parameter 
space of the then available proper motion measurement of 
~ 141kms~^, we did not find any important changes in 
the kick velocity and progenitor constraints related to the 
decrease of the minimum admissible helium star mass. 
With the new proper motion upper limit of 30 kms~^, the 
adopted minimum helium star mass has a much stronger 
effect on the available parameter space. In particular, 
when both the new proper motion constraint and the 
lower minimum helium star mass are imposed, there is no 
longer a lower limit on the magnitude of the kick velocity 
imparted to pulsar B at birth. Consequently a symmetric 
SN explosion becomes a viable formation mechanism for 
pulsar B. Moreover, allowing pre-SN helium star masses 
smaller than 2.1 M0 also removes the lower limit of ^ 
113° on the angle between the kick velocity direction and 
the direction of the helium star's pre-SN orbital velocity, 
so that the kicks no longer have to be directed opposite 
to the orbital motion. 



B. Piran & Shaviv [19, 20] 

More recently, Piran & Shaviv [19, 20] (hereafter PS) 

used a Monte-Carlo method to argue that the proba- 
bility of finding the system close to the Galactic plane 
is highest if pulsar B received only a marginal kick at 
birth of less than 30kms~^, and if its pre-SN progen- 
itor had a very low mass of only about 1.45 M0. In 
their analysis, they considered the post-SN velocity of 
the binary randomly projected on the axis vertical to the 
Galactic plane (the Z-axis introduced in § IV A) and fol- 
lowed this one component of the Galactic motion forward 
in time for 50 Myr (their assumed age for PSRJ0737- 
3039). However, restricting the Galactic motion to the 
Z-component can potentially bias the PDFs, since the 
horizontal velocity component remains entirely uncon- 
strained (effectively, they accept all post-SN peculiar ve- 
locities as equally probable regardless of their horizontal 
component). In particular, if this horizontal component 
is large, it reduces the probability that the system would 
be found at its current X and Y position in a similar 
way as a large vertical component reduces the probabil- 
ity that the system is found at its current Z position. 
The analysis of PS furthermore implicitly assumes that 
the vertical motion is decoupled from the horizontal one. 
Although this is true locally, it is not true globally, since 
the vertical acceleration does not depend just on the ver- 
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FIG. 7: Onc-dimonsional PDFs illustrating some of the dependencies of the derived pulsar B kick velocity and progenitor 
properties on the adopted model assumptions. All plots are for a present-day transverse velocity of lOkms"^. Solid lines 
correspond to uniform radial velocity distributions and dashed lines to Gaussian distributions with velocity dispersions of 
60kms~^ (long-dashed lines) , 130kms~^ (short-dashed line), and 200 km s~^ (dash-dotted line) . For clarity, the PDFs are offset 
from each other by an arbitrary amount. Panels (a)-(c) show the distributions of kick velocity magnitudes Vk, panels (d)-(e) 
the distributions of kick direction cosines cos 9, and panel (f) the distributions of pre-SN helium star masses Mq. 



tical distance from the plane but also on the horizontal 
distance from the Galactic center. Therefore, if the bi- 
nary has traveled through a wide range of Galactocentric 
radii, as indicated by the large range of possible birth 
sites displayed in Fig. 8, its vertical motion is different 
from the one modeled by PS. 

We investigate this more thoroughly by following the 
steps outlined by PS. Firstly, we determine the post- 

SN orbital semi-major axis A and eccentricity e by in- 
tegrating the equations governing the orbital evolution 



due to gravitational wave emission backwards in time, 
for an assumed age of 50Myr. We find A = 1.42 i?0 and 
e = 0.105. 

Secondly, we randomly generate a pulsar B kick veloc- 
ity magnitude Vk, and kick direction angles and 0, as- 
suming a uniform distribution for the magnitude and an 
isotropic distribution for the direction of the kick. Equa- 
tions (4) and (5) and the associated constraints are then 
solved for Mq and Aq, using the values of A and e de- 
rived by reversing the orbital evolution for 50Myr [59]. 
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FIG. 8: Distribution of possible birth sites of PSR J0737-3039 
in the Galactic plane for a transverse velocity component of 
30kms~^, a minimum helium star mass of 1.25 Mq, an age 
range of 30-70 Myr, and a Gaussian radial velocity distribu- 
tion with a velocity dispersion of 130kms~^. The solid circle 
indicates the system's present position. 



If a solution for Mq and Aq exists, the magnitude of the 
post-SN pecuHar velocity imparted to the binary's center 
of mass is calculated by means of Eq. (6). 

Thirdly, we approximate the motion of the system in 
the Galaxy by a local vertical oscillation with an am- 
plitude drawn from a Gaussian distribution with a dis- 
persion of 50 pc. At the time of the SN explosion, the 
system is assumed to be at a uniformly distributed ran- 
dom phase of the oscillation. The vertical motion of the 
system in the Galaxy is then followed forward in time 
for 50 Myr using the Galactic potential of Paczynski [45] . 
For this purpose, the randomly generated starting point 
has to be supplemented with the Z-component of the 
calculated post-SN peculiar velocity. Since the kick di- 
rection angles 9 and (p are defined with respect to the 
pre-SN orbital plane, which has an unknown direction in 
space, the Z-component of the post-SN peculiar velocity 
is determined by assigning an isotropically distributed 
random direction to the calculated post-SN peculiar ve- 
locity. This new direction is generated with respect to the 
Galactic frame of reference OXYZ introduced in § IV A, 
and therefore straightforwardly leads to the determina- 
tion of the Z-component of the post-SN peculiar velocity. 

If the generated system satisfies all SN orbital dynam- 
ics constraints and if, after following the motion in the 
Galaxy for 50 Myr, the system ends up within 50 pc of 
the Galactic mid-plane, it is considered a possible pro- 
genitor of PSRJ0737-3039 [60]. The procedure is then 
repeated until 10* such systems are found. The outcome 
of the simulation is a 2-D PDF of Mq and Vk shown 
in panel (a) of Fig. 9. The PDF is normalized so that 



its integral over the entire admissible parameter space 
is equal to one. This is different from PS who normalize 
their PDF to the peak probability density [61]. The over- 
all features of this PDF are to be compared to the upper 
right panel of Fig. 1 in Piran & Shaviv [19]; it is evident 
that the two PDFs are in rather good agreement. Differ- 
ences in the details can be partly attributed to the PDF 
normalization adopted by PS (since it requires that the 
exact PDF peak is "captured" by the Monte Carlo pro- 
cess), which depends on the accuracy of the Monte-Carlo 
simulation and the size of the bins used to construct the 
PDF. 

In order to assess the effect of approximating the space 
motion of PSR J0737-3039 as a vertical oscillation cen- 
tered on the system's current position, we repeated the 
above calculation for a system completely identical to 
the double pulsar, but located at a different position in 
the Galactic plane. In particular, we positioned the sys- 
tem closer to the Galactic center a.t X = —2.5 kpc and 
Y = —3.5 kpc (this corresponds to the most likely birth 
site found in § VI for a present-day transverse velocity of 
30kms~^, an age of 30-70 Myr, a Gaussian radial veloc- 
ity distribution with a velocity dispersion of 130 kms^^, 
and a minimum pre-SN helium star mass of 2.1 Mq). 
Somewhat surprisingly, we found the effect on the PDFs 
to be rather minimal. As a second test, we therefore con- 
sidered a second system identical to the double pulsar, 
but located at a height of 790 ± 50 pc above the Galactic 
plane (corresponding to the distance of PSRB1534-M2 to 
the Galactic plane). Again, we found minimal changes in 
the PDF. We conclude that the position of the system in 
the Galaxy therefore seems to be relatively unimportant 
in the derivation of the {Vk, Afo)-PDF. 

This conclusion is further strengthened by the results 
in panel (b) of Fig. 9: it shows the 2-D PDF obtained 
without any constraint on |Z|. It is evident that the PDF 
does not depend in any significant way on the position of 
PSR J0737-3039 close to the Galactic plane, in contrast 
to what is claimed by PS (the main difference between 
panels (a) and (b) is that the "peak" contour is somewhat 
more spread out). Indeed when tracking the progenitor 
systems as they are affected by the different constraints 
in the simulations, we find that the \Z\ < 50 pc constraint 
eliminates only 13% of the systems that satisfy all other 
imposed constraints. It turns out that the orbital dy- 
namics of asymmetric SN explosions and associated con- 
straints are much more important than the position of 
the system in the Galaxy in determining the (Vfe,Mo)- 
PDF. 

For comparison, we also constructed 2-D kick velocity 
and progenitor mass distributions by marginalizing the 
5-D PDFs derived in § VI, for model assumptions similar 
to those adopted by Piran & Shaviv [20] . Specifically, we 
considered a present-day transverse velocity of 30 km s~^ , 
a minimum helium star mass of 1.25 Mq, kinematic ages 
of 49-51 Myr, and a Gaussian radial velocity distribution 
with a dispersion of 130 km s"-'^. We also removed the 
limits on pulsar A's spin-orbit misalignment from the list 
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FIG. 9: Contour plots for the (14, Mo)-PDF. Panels (a) and (b) correspond to the PDFs obtained using the method outlined by 
PS. The PDF shown in Panel (a) only includes systems with \Z\ < 50 pc, while the PDF shown in panel (b) does not incorporate 
any constraints on the position of the system in the Galaxy. Panels (c) and (d) correspond to the PDFs obtained using the 
method outlined in this paper. The PDF shown in panel (c) corresponds to a present-day transverse velocity of 30kms~^, a 
minimum pre-SN helium star mass of 1.25 A'Iq, kinematic ages of 49-51 Myr, and a Gaussian radial velocity distribution with a 
dispersion of 130kms~^. The PDF shown in panel (d) adopts the same assumptions, except that the minimum pre-SN helium 
star mass is assumed to be 2.1 Mq. 



of imposed constraints as it was not incorporated in the 
analysis of PS. The resulting PDF is shown in panel (c) 
of Fig. 9. We note that, although the analysis by PS does 
not incorporate the radial velocity, it is an indispensable 
part of our method which cannot be omitted from the 
analysis. A comparison between our results and those de- 
rived by PS with exactly the same assumptions is there- 
fore not possible. Despite this, the differences between 
the PDFs displayed in panels (a) and (c) are quite mini- 
mal. The same holds true if the Gaussian radial velocity 
distribution with a dispersion of 130 km s^^ is replaced 
by any of the other radial- velocity distributions consid- 
ered. In view of the non-negligible differences between 
the two approaches, this agreement is all but trivial and 
thus strengthens our confidence in the validity of the de- 
rived results. The effect of the minimum helium star mass 
on the PDF is illustrated in panel (d), where the mini- 
mum pre-SN helium star mass was assumed to be 2.1 Mq. 
The main difference with the PDF shown in panel (c) is 
the minimum kick velocity of ~ 60kms~^. The most 
likely kick velocity is of the order of 60-100 km s~^ and 
the most likely progenitor mass is 2.1 — 2.4 Mq. 

Thus, we conclude that (i) if one allows a priori for low- 
mass progenitors of pulsar B, and (ii) if the marginalized 



2-D PDF is considered and not the underlying 5-D PDF, 
then the preference for low progenitor masses and low 
kicks is rather generic and it is not related to the sys- 
tem's position close to the plane or the small upper limit 
on the proper motion in any crucial way (c.f. PS). This 
applies to both the analysis presented in this paper and 
the analysis by PS. Examining the role of the prior as- 
sumptions and considering the full 5-D PDF of 14, cos 6, 
4>, Mq, and Aq are therefore essential in establishing an 
unbiased view of the possible formation mechanism and 
progenitor of PSR J0737-3939B. 



Finally, we note that, contrary to what is claimed by 
PS, small tilt angles between pulsar A's spin angular mo- 
mentum and the post-SN orbital angular momentum do 
not necessarily imply that a small SN kick was imparted 
to pulsar B at birth. This can, e.g., be seen from Figs. 
11 and 12 in Paper II. In addition, Eq. (7) shows that the 
tilt depends on the direction as well as on the magnitude 
of the natal kick velocity, and that in the particular case 
where the kicks are confined to the pre-SN orbital plane, 
the tilt angle is equal to zero regardless of the magnitude 
of the kick velocity. 
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VIII. SUMMARY AND CONCLUDING 
REMARKS 

In this study, wc have considered the kinematic and 
evolutionary history of the double pulsar PSRJ0737- 
3039, taking into account the most up-to-date obser- 
vational constraints. Unlike our earlier work on this 
topic, we focus here on the derivation of the full multi- 
dimensional probability distribution functions for the 
magnitude and direction of the kick velocity imparted 
to pulsar B at birth, the mass of pulsar B's pre-SN he- 
lium star progenitor, and the pre-SN orbital separation. 
The consideration of the multi-dimensional PDF has the 
distinct advantage that it is free of projection effects and 
that it accounts for all possible correlations between the 
derived parameters; these correlations and projection ef- 
fects turn out to be very important when addressing the 
question of what are the most likely properties of the dou- 
ble pulsar's progenitor. We also examine the dependence 
of the PDFs on the prior assumptions for the magnitude 
of the transverse systemic velocity component, the mini- 
mum pre-SN helium star mass required to form a NS, the 
age of the system, and the present-day radial velocity. 
Since the latter cannot be determined from observation, 
we construct and explore theoretical radial-velocity dis- 
tributions by means of population synthesis calculations 
for coalescing DNSs. 

One of the main results of our analysis is that marginal- 
izing the full fivc-dimensional progenitor PDF to 1-D or 
2-D distributions for the pulsar B kick velocity and pro- 
genitor mass has a major effect on the determination of 
their most likely values. When the full multi-dimensional 
PDF is examined, it is clear that although some sets of 
prior assumptions indeed favor low kick velocities and 
low progenitor masses, as claimed by Piran & Shaviv 
[19], the majority of the models favor kick velocities of 
50-180 kms~^ and progenitor masses of 1A5-2.75Mq. 

In particular, if the transverse systemic velocity is as- 
sumed to be 30km/s and helium stars less massive than 
2.1 Mq are assumed to be viable NS progenitors, the 
most likely pulsar B progenitor mass is 1.45 regard- 
less of any of the other prior assumptions. If, on the 
other hand, the transverse systemic velocity is assumed 
to be lOkm/s while keeping the assumption that helium 
stars less massive than 2.1 Mq are viable NS progenitors, 
the most likely progenitor mass can vary from 1.35 M^. 
to 2.65 Mq, depending on the assumed systemic ago and 
radial velocity distribution. Hence, whether progenitor 
masses greater than 2.1 Mq are statistically likely or un- 
likely depends strongly on the adopted prior assumptions 
(cf. [19, 20]). Most likely progenitors with Mq < 2.1 Mq 
can furthermore also be associated with kick velocities of 
up to lOOkms"^ (see, e.g., the case of 14 = 30kms~^ 
and Mq > 1.25 Mq in Figs. 4 and 6), while kick veloci- 
ties of less than 60 km s~^ are only possible for progenitor 
masses below 2.1 Mq. 

We also find that the proximity of PSR J0737-3039 to 
the Galactic plane and the small proper motion do not 



pose stringent constraints on the kick velocity and pro- 
genitor mass of pulsar B. Instead, the constraints are pre- 
dominantly determined by the orbital dynamics of asym- 
metric SN explosions. This is in contrast to the work of 
Piran & Shaviv [19, 20] who emphasize that the prox- 
imity of the double pulsar to the Galactic plane implies 
that pulsar B most likely received only a small kick at 
birth and that its progenitor most likely had mass of 

~ 1.45 Mq. 

Hence, based on the currently available observational 
constraints, a wide range of progenitor and kick velocity 
properties are favored for PSR J0737-3039B. In particu- 
lar, the constraints are compatible with a conventional 
hydro dynamical or neutrino-driven SN explosion from a 
helium star more massive than 2.1 Mq [10, 41], as well as 
an electron-capture SN from a helium star less massive 
than 2.5 Mq [50, 51]. Podsiadlowski et al. [3] have spec- 
ulated that the electron-capture SN mechanism may be 
typical for close binaries and that it may be accompanied 
by much smaller kicks than hydrodynamical or neutrino- 
driven SN explosions. Consequently, if pulsar B is formed 
through an electron capture SN and if electron capture 
SNe are accompanied by small kicks, the mass of pul- 
sar B's pre-SN helium star progenitor must be smaller 
than 2.1 Mq (otherwise the kick is always larger than 
~ 60kms~^). Since it is unlikely that future observa- 
tions will lead to new constraints on the smallest possi- 
ble pulsar B progenitor mass, further insight to the for- 
mation mechanism of pulsar B should be sought in core- 
collapse simulations of low-mass (^2.1 Mq) helium stars 
and population synthesis studies of PSR J0737-3039-type 
binaries and their progenitors. 
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APPENDIX A: COMMENTS BY PIRAN &; 
SHAVIV ON THIS WORK 

Piran & Shaviv [21] recently posted a Comment on the 
astro-ph preprint server " refuting" the work presented in 
this paper. Their argumentation, however, rests on se- 
vere misstatements of our analysis and neglects some of 
the major conclusions drawn in the previous sections. 
While we have already addressed their criticism at ap- 
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propriate places throughout this paper, we here wish to 
summarize our counterarguments to ehminate any fur- 
ther misunderstandings: 

(1) Our results are not "in spite of the observations" 
as claimed by Piran & Shaviv [21]. They are based on 
all currently available observational constraints, most of 
which are in fact initial conditions in our calculation. 

(2) While we do conclude that most of our models yield 
kick velocities of 50-180 km s~^, we do not claim that the 
models favoring lower kick velocities are less likely. We 
conclude even less that the most likely pulsar B progen- 
itor mass is large. It can, e.g., be seen from Table IV 
that if the transverse systemic velocity is assumed to be 
30 km/s and the minimum pulsar B progenitor mass is as- 
sumed to be 1.25 A/q, the most likely progenitor mass is 
1.45 Mq regardless of any of the other prior assumptions. 
If, on the other hand, the transverse systemic velocity is 
assumed to be 10 km/s while the minimum progenitor 
mass is kept at 1.25M©, the most likely progenitor mass 
can vary from 1.35 Mq to 2.65 Mq, depending on the as- 
sumed age range and radial velocity distribution. Our 
main conclusion is thus that the prior assumptions play 
an important role in determining the most likely kick ve- 
locity and progenitor properties. 

(3) We do not assume that the system is moving with 
a very large velocity almost exactly towards. Instead, 
we adopt theoretical radial-velocity distributions based 
on realistic binary evolution models for coalescing DNSs 
that reside close to the Galactic plane and have small 
transverse systemic velocity components, consistent with 
the current observational constraints on the kinematics 
of the double pulsar. Because of the latter constraints, 
the distributions implictly incorporate the geometrical 
probability argument put forward in Eq. (2) of [21]. We 
furthermore also consider four different distribution func- 
tions (one flat distribution and three Gaussians of vary- 
ing width) in order to assess the robustness of our results 



(see, e.g.. Table IV). 

(4) The differences between our conclusions and those 
of Piran & Shaviv [19] rest on the consideration of the 
full 5-D PDF instead of a marginalized 2-D one, not on 
the modeling of the motion of the system in Galaxy. The 
reason for this is that the orbital parameters and the 
dynamics of asymmetric SN explosions turn out to be 
much more constraining than the system's position and 
motion in the Galaxy (contrary to what is claimed by 
Piran & Shaviv [19]). Piran & Shaviv [21] also completely 
miss that when we adopt similar assumptions as they do, 
and marginalize our 5-D PDFs for Vk, cos 9, (f), Mq, and 
Aq to 2-D PDFs for Mq and Vk, we find results that 
are in excellent agreement with theirs, irrespective of the 
treatment of the unknown radial velocity. 

(5) The population synthesis calculations used to con- 
struct the theoretical radial-velocity distributions incor- 
porate the possibility that mass transfer takes place right 
before the formation of pulsar B. The models therefore al- 
low the pre-SN masses of helium stars in DNS progenitors 
to be as low as the 1.45 Mq pulsar B progenitor mass fa- 
vored by Piran & Shaviv [19]. The derived radial- velocity 
distributions are thus by no moans biased towards high 
progenitor masses as claimed by Piran & Shaviv [21]. 

(6) Based on the preceding arguments, the conclusion 
of Piran & Shaviv [21] that "we assume the result we 
obtain" clearly rests on severe misinterpretations of the 
calculations underlying our analysis. 

We want to conclude this discussion by stressing that 
the goal of our paper is not to refute the progenitor and 
kick constraints derived by Piran & Shaviv [19]. We do 
however want to make it clear that there are many addi- 
tional factors that were not considered by these authors, 
and that these may affect the determination of the most 
likely progenitor and kick velocity properties. 
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orbital separation, the magnitude and direction of pul- 
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motion preceding the second supernova explosion is most 
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our radial velocity distributions with and without this 
early motion and we can confirm that its effect is in- 
deed negligible. This is not so much due to the somewhat 
smaller center-of-mass velocities after the first explosion 
(compared to those after the second explosion because of 
the wide orbits involved); instead it is mostly due to the 
very short time between between the two events (from 
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sar B's pre-SN helium star progenitor and the small val- 
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results presented in this paper. 

[58] Confidence limits for cos6', (f>, and Aq were also calcu- 
lated but these did not contribute any significant new 
information. 

[59] Note that PS consider a range of eccentricities 0.088 < 
e < 0.14 obtained by reversing the past orbital evolution 
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acteristic age of pulsar A). Although, this is inconsistent 
with the age of 50 Myr adopted for the Galactic motion 
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For these reasons, we modify this step in PS's analysis 
and consistently use the post-SN orbital parameters as- 
sociated with a fixed age of 50 Myr. 
[60] Contrary to PS, we do not impose any constraint on the 
present-day transverse systemic velocity, since there is no 
unambiguous way to impose the constraint without mak- 
ing additional assumptions about the system's present- 
day 3-D space velocity. In particular, limiting the motion 
of the system in the Galaxy to a local vertical oscillation 
restricts the knowledge of the velocity in the simulations 



to the component perpendicular to the Galactic plane. 
Since this component does not lay in the plane perpen- 
dicular to the linc-of-sight to the double pulsar, imposing 
the proper motion constraint requires additional assump- 
tions on the radial velocity and the direction of the proper 
motion. 

[61] Because of the chosen normalization, the contour levels 
in PS measure the decrease of the PDF with respect to 
its maximum value. They therefore do not correspond to 
confidence levels. 



